#  Project: War Time Military Service Can Affect Partisan Preferences
#  Last updated: November 8,2022
#  Purpose: Create main Figures and Tables for Paper
#  Outputs: data/ines_77_2019.csv"
#  Machine: Chagai's macbook pro

#Load Packages------
library("tidyverse")
library("ggplot2")
library("stargazer")
library("xtable")
library("texreg")
library("readstata13")
library("estimatr")
library("naniar")
library("lfe")
library("multiwayvcov")
library("lmtest")
library("rdrobust")
library("rddensity")
library("sjmisc")
library("rdd")



#Read Data------
ines_master <- read.csv("data_for_analyses/ines_77_2019.csv") %>% 
  filter(.,
         aliya_after_73 !=1,
         age_in_73 < 74)
ines_master_olim <- read.csv("data_for_analyses/ines_77_2019.csv")# data including post 1973 immigrants


#Table 1-----

# Create Main Table

labor <- rdrobust(ines_master$vote_labor, 
                  ines_master$age_in_73,
                  c = 18, p = 1)
summary(labor)


likud <- rdrobust(ines_master$vote_likud, 
                  ines_master$age_in_73,
                  c = 18, p = 1)
summary(likud)


undecided <- rdrobust(ines_master$undecide, 
                      ines_master$age_in_73,
                      c = 18, p = 1)
summary(undecided)


turnout <- rdrobust(ines_master$no_vote, 
                    ines_master$age_in_73,
                    c = 18, p = 1)
summary(turnout)


# Create Data frame for Table
main_result <- rbind(
  cbind(round(labor$coef,3), round(labor$se,3), 
        rbind(round(labor$bws[1,1]),round(labor$bws[1,1]),round(labor$bws[1,1])), 
        rbind(labor$p,labor$p,labor$p), sum(labor$N), round(labor$ci,3)),
  cbind(round(likud$coef,3), round(likud$se,3), 
        rbind(round(likud$bws[1,1]),round(likud$bws[1,1]),round(likud$bws[1,1])), 
        rbind(likud$p,likud$p,likud$p), sum(likud$N), round(likud$ci,3)),
  cbind(round(undecided$coef,3), round(undecided$se,3), 
        rbind(round(undecided$bws[1,1]),round(undecided$bws[1,1]),round(undecided$bws[1,1])), 
        rbind(undecided$p,undecided$p,undecided$p), sum(undecided$N), round(undecided$ci,3)),
  cbind(round(turnout$coef,3), round(turnout$se,3), 
        rbind(round(turnout$bws[1,1]),round(turnout$bws[1,1]),round(turnout$bws[1,1])), 
        rbind(turnout$p,turnout$p,turnout$p), sum(turnout$N), round(turnout$ci,3)))

main_result_df <- as.data.frame(main_result) %>% 
  dplyr :: rename(.,
                  "Estimate" = "Coeff",
                  "SE" = "Std. Err.",
                  "Bandwidth" = "V3",
                  "Poly" = "V4",
                  "Obs." = "V5") %>% 
  slice(3,6,9,12) %>% 
  rownames_to_column() %>% 
  mutate(.,
         rowname = case_when(
           rowname == "Robust" ~ "Labor",
           rowname == "Robust.1" ~ "Likud",
           rowname == "Robust.2" ~ "Undecided",
           rowname == "Robust.3" ~ "No Vote"
         ),
         rename(.,
                Outcome = rowname)) %>% 
  select(.,
         -Poly, -"CI Lower", -"CI Upper", -Outcome)

main_result_df<- rotate_df(main_result_df, cn =T)

main_result_df[2,] <- paste0("(", format(unlist(main_result_df[2,])),")")  
rownames(main_result_df)[rownames(main_result_df) == "Estimate"] <- ""
rownames(main_result_df)[rownames(main_result_df) == "SE"] <- " "
stargazer(main_result_df,
          summary = FALSE, rownames = T,
          title = "RD Estimates -- The Effects of Participation in 1973 War",
          label = "tab:rd_1",
          column.sep.width = "15pt",
          style = "ajps",
          font.size = "small",
          header=FALSE, 
          digit.separator	= ",",
          digits = 3,
          notes = c("Regression discontinuity models with MSE optimal, bandwidths",
                    "and a triangular kernal. Robust standard errors in parentheses."),
          notes.append = F,
          notes.align = "l")



# Figure 1------
rdplot(ines_master$vote_labor, 
       ines_master$age_in_73,  
       c = 18, p = 1,
       x.label = "Age in 1973", 
       y.label = "Vote Labor",
       title = "")

#Figure 2------
labor_all <- rdrobust(ines_master$vote_labor, 
                      ines_master$age_in_73,
                      c = 18, p = 1)
summary(labor_all)

labor_post77 <- rdrobust(ines_master$vote_labor[ines_master$election_year>1977], 
                         ines_master$age_in_73[ines_master$election_year>1977],
                         c = 18, p = 1)
summary(labor_post77)

labor_post81 <- rdrobust(ines_master$vote_labor[ines_master$election_year>1982], 
                         ines_master$age_in_73[ines_master$election_year>1982],
                         c = 18, p = 1)
summary(labor_post81)

labor_post84 <- rdrobust(ines_master$vote_labor[ines_master$election_year>1984], 
                         ines_master$age_in_73[ines_master$election_year>1984],
                         c = 18, p = 1)
summary(labor_post84)

labor_post88 <- rdrobust(ines_master$vote_labor[ines_master$election_year>1989], 
                         ines_master$age_in_73[ines_master$election_year>1989],
                         c = 18, p = 1)
summary(labor_post88)

labor_post92 <- rdrobust(ines_master$vote_labor[ines_master$election_year>1992], 
                         ines_master$age_in_73[ines_master$election_year>1992],
                         c = 18, p = 1)
summary(labor_post92)


labor_post96 <- rdrobust(ines_master$vote_labor[ines_master$election_year>1996], 
                         ines_master$age_in_73[ines_master$election_year>1996],
                         c = 18, p = 1)
summary(labor_post96)

labor_post99 <- rdrobust(ines_master$vote_labor[ines_master$election_year>1999], 
                         ines_master$age_in_73[ines_master$election_year>1999],
                         c = 18, p = 1)
summary(labor_post99)

labor_post01 <- rdrobust(ines_master$vote_labor[ines_master$election_year>2001], 
                         ines_master$age_in_73[ines_master$election_year>2001],
                         c = 18, p = 1)
summary(labor_post01)


labor_post03 <- rdrobust(ines_master$vote_labor[ines_master$election_year>2003], 
                         ines_master$age_in_73[ines_master$election_year>2003],
                         c = 18, p = 1)
summary(labor_post03)

labor_post06 <- rdrobust(ines_master$vote_labor[ines_master$election_year>2006], 
                         ines_master$age_in_73[ines_master$election_year>2006],
                         c = 18, p = 1)
summary(labor_post06)


labor_post09 <- rdrobust(ines_master$vote_labor[ines_master$election_year>2009], 
                         ines_master$age_in_73[ines_master$election_year>2009],
                         c = 18, p = 1)
summary(labor_post09)

labor_post13 <- rdrobust(ines_master$vote_labor[ines_master$election_year>2013], 
                         ines_master$age_in_73[ines_master$election_year>2013],
                         c = 18, p = 1)
summary(labor_post13)


labor_post15 <- rdrobust(ines_master$vote_labor[ines_master$election_year>2015], 
                         ines_master$age_in_73[ines_master$election_year>2015],
                         c = 18, p = 1)
summary(labor_post15)




long_term <- rbind(cbind(labor_all$coef, labor_all$ci, sum(labor_all$N), "1977-2019", labor_all$pv),
                   cbind(labor_post77$coef, labor_post77$ci, sum(labor_post77$N), "1981-2019", labor_post77$pv),
                   cbind(labor_post81$coef, labor_post81$ci, sum(labor_post81$N), "1984-2019", labor_post81$pv),
                   cbind(labor_post84$coef, labor_post84$ci, sum(labor_post84$N), "1988-2019", labor_post84$pv),
                   cbind(labor_post88$coef, labor_post88$ci, sum(labor_post88$N), "1992-2019", labor_post88$pv),
                   cbind(labor_post92$coef, labor_post92$ci, sum(labor_post92$N), "1996-2019", labor_post92$pv),
                   cbind(labor_post96$coef, labor_post96$ci, sum(labor_post96$N), "1999-2019", labor_post96$pv),
                   cbind(labor_post99$coef, labor_post99$ci, sum(labor_post99$N), "2001-2019", labor_post99$pv),
                   cbind(labor_post01$coef, labor_post01$ci, sum(labor_post01$N), "2003-2019", labor_post01$pv),
                   cbind(labor_post03$coef, labor_post03$ci, sum(labor_post03$N), "2006-2019", labor_post03$pv),
                   cbind(labor_post06$coef, labor_post06$ci, sum(labor_post06$N), "2009-2019", labor_post06$pv),
                   cbind(labor_post09$coef, labor_post09$ci, sum(labor_post09$N), "2013-2019", labor_post09$pv))

long_term <- as.data.frame(long_term) %>% 
  slice(3,6,9,12,15,18,21,24,27,30,33,36,39,42) %>% 
  rename(.,
         low = "CI Lower",
         high = "CI Upper") %>% 
  mutate(.,
         low = as.numeric(low),
         high = as.numeric(high), 
         Years = as.character(V5), 
         Coeff = as.numeric(Coeff))

# Plot long term
ggplot(long_term, aes(x = V5, y = Coeff)) +
  geom_hline(yintercept = 0, color = "gray50", linetype = 2, size = 0.2) +
  geom_pointrange(aes(ymin = low, ymax = high), color = "dodgerblue4") +
  ylim(-.35,.14)+
  scale_fill_manual(values =  c("dodgerblue2")) +
  annotate("text", x = 1:12, y = -0.337, 
           label = c("N = 16,962","15,596", "12,153", "10,995","9,796",
                     "8,723","7,664", "6,728","5,798", "4,988",
                     "4,081", "3,309"),#, "2,141", "1,087"),
           color = "gray30", size = 2.8) +
  annotate("text", x = 1:12, y = 0.12, 
           label = c("p = 0.03","0.03", "0.05", "0.09","0.1",
                     "0.14", "0.4","0.18", "0.19",
                     "0.09", "0.15", "0.15"),#, "0.13", "0.35"),        
           color = "gray30", size = 2.8) +
  labs(x = "Election Years Included in RDD",
       y = "RD Effect on Support for Labor") +
  theme(panel.grid.major = element_blank(), 
        axis.text.x = element_text(size = 8, angle = 90),
        plot.caption = element_text(size = 10, family = "Times",hjust = -.02),
        panel.grid.minor = element_blank(),
        panel.background = element_blank(), 
        axis.line = element_line(colour = "black"))


#Table 2-----
# Ideology
ideo <- rdrobust(ines_master$ideology_7, 
                 ines_master$age_in_73,
                 c = 18, p = 1)
summary(ideo)

right <- rdrobust(ines_master$right_wing, 
                  ines_master$age_in_73,
                  c = 18, p = 1)
summary(right)

left <- rdrobust(ines_master$left_wing, 
                 ines_master$age_in_73,
                 c = 18, p = 1)
summary(left)


# Personal effects

spending <- rdrobust(ines_master$spending, 
                     ines_master$age_in_73,
                     c = 18, p = 1)
summary(spending)

rooms <- rdrobust(ines_master$rooms_hh, 
                  ines_master$age_in_73,
                  c = 18, p = 1)
summary(rooms)

ppl <- rdrobust(ines_master$ppl_hh, 
                ines_master$age_in_73,
                c = 18, p = 1)
summary(ppl)

# Evaluating the party

lab_sec <- rdrobust(ines_master$labor_sec, 
                    ines_master$age_in_73,
                    c = 18, p = 1)
summary(lab_sec)


lab_econ <- rdrobust(ines_master$labor_econ, 
                     ines_master$age_in_73,
                     c = 18, p = 1)
summary(lab_econ)

wop <- rdrobust(ines_master$war_over_peace, 
                ines_master$age_in_73,
                c = 18, p = 1)
summary(wop)

# Create Data frame for Table
addition_out <- rbind(
  cbind(round(lab_sec$coef,3), round(lab_sec$se,3), 
        rbind(round(lab_sec$bws[1,1]),round(lab_sec$bws[1,1]),round(lab_sec$bws[1,1])), 
        rbind(lab_sec$p,lab_sec$p,lab_sec$p), sum(lab_sec$N), round(lab_sec$ci,3)),
  cbind(round(wop$coef,3), round(wop$se,3), 
        rbind(round(wop$bws[1,1]),round(wop$bws[1,1]),round(wop$bws[1,1])), 
        rbind(wop$p,wop$p,wop$p), sum(wop$N), round(wop$ci,3)),
  cbind(round(ideo$coef,3), round(ideo$se,3), 
        rbind(round(ideo$bws[1,1]),round(ideo$bws[1,1]),round(ideo$bws[1,1])), 
        rbind(ideo$p,ideo$p,ideo$p), sum(ideo$N), round(ideo$ci,3)),
  cbind(round(right$coef,3), round(right$se,3), 
        rbind(round(right$bws[1,1]),round(right$bws[1,1]),round(right$bws[1,1])), 
        rbind(right$p,right$p,right$p), sum(right$N), round(right$ci,3)),
  cbind(round(left$coef,3), round(left$se,3), 
        rbind(round(left$bws[1,1]),round(left$bws[1,1]),round(left$bws[1,1])), 
        rbind(left$p,left$p,left$p), sum(left$N), round(left$ci,3)),
  cbind(round(spending$coef,3), round(spending$se,3), 
        rbind(round(spending$bws[1,1]),round(spending$bws[1,1]),round(spending$bws[1,1])), 
        rbind(spending$p,spending$p,spending$p), sum(spending$N), round(spending$ci,3)),
  cbind(round(rooms$coef,3), round(rooms$se,3), 
        rbind(round(rooms$bws[1,1]),round(rooms$bws[1,1]),round(rooms$bws[1,1])), 
        rbind(rooms$p,rooms$p,rooms$p), sum(rooms$N), round(rooms$ci,3)),
  cbind(round(ppl$coef,3), round(ppl$se,3), 
        rbind(round(ppl$bws[1,1]),round(ppl$bws[1,1]),round(ppl$bws[1,1])), 
        rbind(ppl$p,ppl$p,ppl$p), sum(ppl$N), round(ppl$ci,3)))

addition_out_df <- as.data.frame(addition_out) %>% 
  dplyr :: rename(.,
                  "Estimate" = "Coeff",
                  "SE" = "Std. Err.",
                  "Band" = "V3",
                  "Poly" = "V4",
                  "Obs." = "V5") %>% 
  slice(3,6,9,12,15,18,21,24,27) %>% 
  rownames_to_column() %>% 
  mutate(.,
         rowname = case_when(
           rowname == "Robust" ~ "Labor Security",
           rowname == "Robust.1" ~ "Prepare War",
           rowname == "Robust.2" ~ "Ideology",
           rowname == "Robust.3" ~ "Right",
           rowname == "Robust.4" ~ "Left",
           rowname == "Robust.5" ~ "Spending",
           rowname == "Robust.6" ~ "HH Rooms",
           rowname == "Robust.7" ~ "HH Members"),
         rename(.,
                Outcome = rowname)) %>% 
  select(.,
         -Poly, -"CI Lower", -"CI Upper", -Outcome)

addition_out_df<- rotate_df(addition_out_df, cn =T)

addition_out_df[2,] <- paste0("(", format(unlist(addition_out_df[2,])),")")  
rownames(addition_out_df)[rownames(addition_out_df) == "Estimate"] <- ""
rownames(addition_out_df)[rownames(addition_out_df) == "SE"] <- " "

stargazer(addition_out_df,
          summary = FALSE, rownames = T,
          title = "RD Estimates -- Potential Mechanisms",
          #  dep.var.caption="Males ................... Females",
          label = "tab:rd_add_out",
          column.sep.width = "-4pt", 
          font.size = "small",
          style = "ajps",
          header=FALSE, 
          digit.separator	= ",",
          digits = 3,
          notes = c("Regression discontinuity models with MSE optimal bandwidths, and a triangular kernal. Robust",
                    "standard errors in parentheses."),
          notes.append = F,
          notes.align = "l")


# Figure 3-----
rdplot(ines_master$labor_sec, ines_master$age_in_73,  
       c = 18, p = 1,
       x.label = "Age in 1973", 
       y.label = "Labor is Most Suited to Deal with Security Affairs",
       title = "")


# Figure 4-----
rdplot(ines_master$war_over_peace, ines_master$age_in_73,  
       c = 18, p = 1,
       x.label = "Age in 1973", 
       y.label = "Prepare for War is Best Way to Secure Israel",
       title = "")


#Table 3-----
# Create 67 Table
labor67 <- rdrobust(ines_master$vote_labor, 
                    ines_master$age_in_67,
                    c = 18, p = 1)
summary(labor67)


likud67 <- rdrobust(ines_master$vote_likud, 
                    ines_master$age_in_67,
                    c = 18, p = 1)
summary(likud67)


undecided67 <- rdrobust(ines_master$undecide, 
                        ines_master$age_in_67,
                        c = 18, p = 1)
summary(undecided67)


turnout67 <- rdrobust(ines_master$no_vote, 
                      ines_master$age_in_67,
                      c = 18, p = 1)
summary(turnout67)


# Create Data frame for Table
main_result67 <- rbind(
  cbind(round(labor67$coef,3), round(labor67$se,3), 
        rbind(round(labor67$bws[1,1]),round(labor67$bws[1,1]),round(labor67$bws[1,1])), 
        rbind(labor67$p,labor67$p,labor67$p), sum(labor67$N), round(labor67$ci,3)),
  cbind(round(likud67$coef,3), round(likud67$se,3), 
        rbind(round(likud67$bws[1,1]),round(likud67$bws[1,1]),round(likud67$bws[1,1])), 
        rbind(likud67$p,likud67$p,likud67$p), sum(likud67$N), round(likud67$ci,3)),
  cbind(round(undecided67$coef,3), round(undecided67$se,3), 
        rbind(round(undecided67$bws[1,1]),round(undecided67$bws[1,1]),round(undecided67$bws[1,1])), 
        rbind(undecided67$p,undecided67$p,undecided67$p), sum(undecided67$N), round(undecided67$ci,3)),
  cbind(round(turnout67$coef,3), round(turnout67$se,3), 
        rbind(round(turnout67$bws[1,1]),round(turnout67$bws[1,1]),round(turnout67$bws[1,1])), 
        rbind(turnout67$p,turnout67$p,turnout67$p), sum(turnout67$N), round(turnout67$ci,3)))

main_result67_df <- as.data.frame(main_result67) %>% 
  dplyr :: rename(.,
                  "Estimate" = "Coeff",
                  "SE" = "Std. Err.",
                  "Bandwidth" = "V3",
                  "Poly" = "V4",
                  "Obs." = "V5") %>% 
  slice(3,6,9,12) %>% 
  rownames_to_column() %>% 
  mutate(.,
         rowname = case_when(
           rowname == "Robust" ~ "Labor",
           rowname == "Robust.1" ~ "Likud",
           rowname == "Robust.2" ~ "Undecided",
           rowname == "Robust.3" ~ "No Vote"
         ),
         rename(.,
                Outcome = rowname)) %>% 
  select(.,
         -Poly, -"CI Lower", -"CI Upper", -Outcome)

main_result67_df<- rotate_df(main_result67_df, cn =T)

main_result67_df[2,] <- paste0("(", format(unlist(main_result67_df[2,])),")")  
rownames(main_result67_df)[rownames(main_result67_df) == "Estimate"] <- ""
rownames(main_result67_df)[rownames(main_result67_df) == "SE"] <- " "
stargazer(main_result67_df,
          summary = FALSE, rownames = T,
          title = "RD Estimates -- The Effects of Conscription Eligibility in 1967",
          label = "tab:rd_67",
          column.sep.width = "15pt",
          font.size = "small",
          style = "ajps",
          header=FALSE, 
          digit.separator	= ",",
          digits = 3,
          notes = c("Regression discontinuity models with MSE optimal bandwidths,",
                    "and a triangular kernal. Robust standard errors in parentheses."),
          notes.append = F,
          notes.align = "l")


# Figure 5-----
ines_master1 <- ines_master %>% 
  mutate(.,
         year = case_when(
           election_year == 1981.1 ~ 1981,
           election_year == 1981.2 ~ 1981,
           election_year == 1988.1 ~ 1988
         ),
         year = ifelse(is.na(year), election_year, year))
table(ines_master1$year)


s <- ines_master1 %>% 
  group_by(year)

s_yt <- ines_master1 %>% 
  group_by(year, Sex)

by_year <- summarize(s,
                     labor = mean(vote_labor, na.rm = T),
                     likud = mean(vote_likud, na.rm = T),
                     left = mean(left, na.rm = T),
                     right = mean(right_wing, na.rm = T),
                     center = mean(center, na.rm = T))


by_year_treat <- summarize(s_yt,
                           labor = mean(vote_labor, na.rm = T),
                           likud = mean(vote_likud, na.rm = T),
                           left = mean(left, na.rm = T),
                           right = mean(right_wing, na.rm = T),
                           center = mean(center, na.rm = T))

ggplot() + 
  geom_line(data = by_year, aes(x = year, y = labor), color = "blue") +
  geom_line(data = by_year, aes(x = year, y = likud), color = "red") +
  geom_point(data = by_year, aes(x = year, y = labor), color = "blue") +
  geom_point(data = by_year, aes(x = year, y = likud), color = "red") +
  ylim(0, 0.8)+
  xlab('Year') +
  ylab('% of Respondents Voting for Labor (Blue) and Likud (Red)') +
  theme(panel.grid.major = element_blank(), 
        axis.text.x = element_text(size = 12),
        plot.caption = element_text(size = 10, family = "Times",hjust = -.02),
        panel.grid.minor = element_blank(),
        panel.background = element_blank(), 
        axis.line = element_line(colour = "black"))


